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Multipole matrix elements of Green function of Laplace equation are calculated. The multipole 
matrix elements of Green function in electrostatics describe potential on a sphere which is produced 
by a charge distributed on the surface of a different (possibly overlapping) sphere of the same radius. 
The matrix elements are defined by double convolution of two spherical harmonics with the Green 
function of Laplace equation. The method we use relies on the fact that in the Fourier space the 
double convolution has simple form. Therefore we calculate the multipole matrix from its Fourier 
transform. An important part of our considerations is simplihcation of the three dimensional Fourier 
transformation of general multipole matrix by its rotational symmetry to the one-dimensional Hankel 
transformation. 
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I. INTRODUCTION 

Metallic spheres in a vacuum, minute spherical particles in a fluid, and spherical inclusions in a solid body are 
often considered in statistical physics of dispersive media [ij. Equations which describe those systems are linear and 
have spherical symmetry. It is the case of Laplace equation in dielectrics 0 , Stokes equations in suspensions Q and 
Lame equations describing solid body. The first step in considerations of dispersive media is often to solve a single 
particle problem. In dielectrics it means to find distribution of charge on the surface of a single metallic sphere in 
an external electrostati(^otential. To find the solution of single particle problem it is convenient to take spherical 
symmetry into account To this end one introduces a set of scalar functions on the sphere which is invariant under 
rotations. That is how spherical harmonics Yim (?) enter calculations. 

To pass from considerations of a single particle to the analysis of many particles it demands to answer the following 
question. What is the electrostatic potential produced by a charge distributed in one spherical surface in the area 
occupied by a different sphere. The answer to the above question can be inferred from the following multipole matrix 
element 


= [ (fir f (fir'-d (|r - R| - a) (j)*^ (r - R) G (r - r') -d (|r' - R'| - a) (j)i>m' (r' - R'). 
In the above definition G (r) is the Green function of Laplace equation 0, 


G(r) 


I 

Anr ’ 


( 1 ) 


a is radius of particles, and one-dimensional Dirac delta distribution S (x) is used. Moreover (pim (r) are solid harmonics 
defined by 

Hm (r) = fiYim (f) , 

with spherical harmonics Yim (?) numbered by order I = 0,1,... and azimuthal number m = — Z,..., Z In our 
notation an argument of spherical harmonics is a versor ? {9, (f>) in direction described by angles 0, cj) in spherical 
coordinates. Dirac delta distributions in equation o reduce three-dimensional integrations to integrals over the 
surface of the spheres with radius a centered at positions R and R'. 

We can differentiate between two situations. The first situation corresponds to nonoverlapping spheres, i.e. when 
|R- R' I > 2a. In this case the matrix elements defined by equation m can be inferred from the results in the 
literature Q. They have application e.g. in numerical simulations. The second situation is the case of overlapping 
configurations, |R — R'| < 2a. Even if the particles in the system cannot overlap there may appear a need of calculation 
of overlapping configurations of the multipole matrix elements [G (R, R0]im I'm'- example, it has been recognized 
that overl^ping configurations appear in microscopic explanation of the famous Clauius-Mossotti formula for dielectric 
constant [g- In this context integral /|R^_R^/|< 2 a (R’ RO] ^XP (~*^R0im I'm' multipole numbers I, V has 

been considered. Therefore the multipole matrix elements [G (R, R0]/m Vm’ overlapping configurations play an 
important role in statistical physics considerations of dispersive media. 

In this article we give general expression for the multipole matrix elements [G (R, RO];™ I'm' defined by equation 
dm- The new contribution of the current article is the calculation of [G (R, R')];^ for overlapping configurations 
|R-R' I < 2a. We are going to use the result for overlapping configurations in our statistical physics considerations 
of transport properties of dispersive media in further work. 

The method of calculation of the multipole matrix elements defined by equation (ED is the following. We observe 
that [G (R, ROljm I'm' form of a double convolution of three functions - two solid harmonics and the Green 

function. Therefore we calculate the Fourier transform of the three functions, take their product, and then perform the 
inverse Fourier transform to obtain G (R, R'). An important element of our calculations is to use spherical symmetry 
which allows to reduce the three-dimensional Fourier transform of the multipole matrix to the one-dimensional Hankel 
transform. 


II. MULTIPOLE MATRIX IN THE FOURIER SPACE 

The aim of this article is calculation of integral ED- We can simplify it using homogeneity of the Laplace equation 
which implies that the multipole matrix [G (R, ROJjm I'm' depends on the relative positions R — R^ Therefore from 
now on we consider multipole matrix [G (R)];^ ;/^/ defined by 
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[G(R-R')]z™,zw = [G(R,R')],„.,w 

depending on the relative positions. 

Next, we observe that the multipole matrix given by the formula ([1]) has a form of a double convolution of three 
functions. The three functions are w/'m' (r) = <5 (|r| — a) 4>i'm' (r) /a, conjugated function (r), and the Green 
function G (r). The Fourier transform of a double convolution of the three functions is given by the product of their 
Fourier transforms, therefore 


G(k) 




= (hL(k)G(k)(hjw (k). 


In our calculations we use the following definition of the three-dimensional Fourier transform 


(k) — J d R exp ( zkR.) Glm^Um' (R) : 

with the inverse transformation given by the formula 


Glm.l'm' (R) — 


1 


(2^)^ 


d?k exp (zkR) Gim,i'm> (k). 


In this situation the Fourier transforms of uiim (r) and G (r) are given respectively by 


G(k) = 


fc2’ 


and 


( 2 ) 


djim (k) = 47ra'+^ (-if ji {ka) Yim (k) ■ 

The last expression can be calculated with the use of equation (5.8.3) from Ref. 0 and with the use of orthonormality 
of spherical harmonics. The expression contains spherical Bessel functions (a;) of the order 1. Finally, the Fourier 
transform of multipole matrix is given by the following formula 


G(k) 


Im.V' 


= ^W'+2 Ji ika)3v {ka) 


(k) Yi,m' (k) . 


( 3 ) 


According to the above expression, the Fourier transform of [G (R)];,„ is given by the spherical harmonics and 
the spherical Bessel functions. At this point we face the main difficulty in our calculations of [G (R)]j,^ . The 

inverse Fourier transform of the above formula has to be calculated. To this end we consider rotational symmetry of 
the multipole matrix. 


III. ROTATIONAL SYMMETRY OF A MULTIPOLE MATRIX 


In the definition of the multipole matrix m there appear solid harmonics 4‘im (j) and Green function G(r). 
Transformation of solid harmonics under rotation is given by the formula 


(t>im (D (a,7,/3) r) 


E 


(a,7,/3)] (r) 

L J mm^ 


( 4 ) 


The above expression can be inferred from Ref. 0 from which we adopt notation in this article. There, formula 
(4.1.1) defines three-dimensional rotation matrix D(a,/ 3 , 7 ) characterized by the Euler angles a, ,5, 7 . Moreover, 
^mm> 7 ) denotes the Wigner matrix. Isotropy of the Laplace equation implies that its Green function G (r) is 
invariant under rotation, i.e. G (D (a,^, 7 )r) = G(r). 
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Changing the variables in the integrals in Eq. (HD and the above properties of the solid harmonics and the Green 
function lead us to the following transformation of the multipole matrix elements under rotation 


Glm,Vm' (^7 T) — 


E 

mi ,mj 






G 


Imi 


(R) • 


(5) 


It is worth mentioning here that the above transformation applied for R in z direction, R = -Re^, and for any rotations 
around axis z, thus for f3 = 0 and any other Euler angles, implies that the multipole matrix is diagonal in indexes m, 
i.e. 


Glm^l'm' — ^m,m'Glm,l'm (R^^) : 

where Kronecker delta Sm,m' appears. It is easy to prove similar diagonality in case of the Fourier transform, 

Glm,l'm' — ^m, 7 n'Gljri,l'm (^®z) : ( 6 ) 

which appears as a result of Fourier transformation of equation ([5]) and considerations for proper rotations and wave 
vector k = 


IV. FOURIER TRANSFORM OF MULTIPOLE MATRIX - SIMPLIFICATION BY ROTATIONAL 

SYMMETRY 


The key point of our calculations is simplification of the Fourier transform of a matrix satisfying symmetry property 
given by equation ([5]) . We simplify the Fourier transform for a wave vector along z direction because then the multipole 
matrix is diagonal in indexes m, as is shown by equation ([HD. For the case k = expression (H)) written in spherical 
coordinates reads 

pCiO pTZ p‘ 2 'K 

Gim,i'fn'(ke^) = / dR d6 d(j)R^ sm6exp{-ikRcos9)Gim,i'fn'(R-{R,0,4>)) ■ (7) 

Jo Jo Jo 

We express vector R (R, 9, (j)) by a product of the vector Re^ and the rotation matrix D (a, /3, 7 ) characterized with 
proper Euler angles. The angles can be deduced from formula (4.1.1) in Ref. [^. These angles are (3 = —9, 7 = —</>, 
and any a, e.g. a = 0. Therefore R (R, 9,(j)) = U (0, —9, —cj)) Re^. For this rotation and R =Rez we use symmetry 
property (HJ which leads to the following expression for multipole matrix G 


Gi™.iw(R(R, 0 , <(>))= ^ 

mi ,mj 




m'm\ V ’ ’ 


Glmi^l'm'^ {R^z') • 


( 8 ) 


In the next step we represent exp {—ikR cos 9) from expression ([71) in the form of an infinite series of spherical Bessel 
functions 

00 

exp {—ikRcos9) = ^ (—0*^ (2^i + 1) jh (kR) Rqo ^ (0, —0, ~4 ’), number 

li—O 

which is deduced from formula (5.8.1) and (4.1.26) in the reference Taking into consideration the last two formulae 
in expression © yields 


pOO OO 

Gim,i'm' {k&z) = / dRR^ ^ {‘Jh + 1) jq {kR) Gimi,i'm[ (Rgz) x 

Zi=0 

pTT plTT r / .N 1 * 

Y. d9 d4> sin ( 0 , 9,4>) ( 0 , 9,4>) dIJ. ( 0 , 9,4>) ■ 

, ^0 ^0 L . 


mi ,mj^ 


(9) 


Integration over variables 9, (p in the above formula is performed with the use of formulae (4.6.2), (4.1.12) and (4.2.7) 
from the reference 0. They lead to expression 




dp sin 9Dqq^ 


{0,9,p)Dg\^^ {0,9,h) 


, ( 0 , 6 »,<))) 
m'.m\ V ’ ’ 


* 


= 4^(-l)’"' 


I I' h 

m —m' 0 


I I' li\ 
mi —m[ 0 J ’ 
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which contains Wigner 3-j symbols. Taking into account the above integral in expression (|ni) yields 




(fee.) = 4^ E E E + 1) 

— — mi—— I m'-^ ——V 

nOO 

{-it / dRRti, {kR) {Re,). 

Jo 


I r h 

m —m! 0 


I 


I' h 


mi —m'l 0 


( 10 ) 


In this way the three-dimensional Fourier transform is reduced to the one dimensional Hankel transform @ . 

It is convenient to introduce different representation of matrix G. Namely, instead of Gim,i'm' {Re,) we can consider 
gf {R) defined in the following way 


gt {R) = (2j + 1) ^ (-I)’" J ) Gim,i'm' {Re,). 

m,m' ^ '' 


( 11 ) 


with the inverse transformation 


{Re,) = (-I)’” E (L-L'o) (^) • (12) 

Let us notice that only integer j which satisfy condition \l — I'\ < j < I + I' needs to be considered in the above 
equations. It follows from properties of Wigner 3-j symbols. The multipole matrix elements Gim,i'm' (R) for any R 
are then related to gf {R) by equation 


Glm,l'm' (R {R, 0, 0)) = ^ (-1)™ 

mi 


l+l' 

E 


I v 


J 


m —m mi 


EE.o(O,-0,-<)>)5/,r (R) 


(13) 


which follows from equation ([5]) and relation m- 

Different representations of the multipole matrix Gim Vm' {Re,) in positional space given by equations dill) and 
m can be similarly introduced also in Fourier space, i.e. 


gt {k) = (2j 1) ^ (-I)’” r I {ke,). 

m,m' ^ ' 


(14) 


Gim.i'm'(fce^,) = (-1)™ ^ (ri -to' o)^D'(*)- (1®) 

By equations (ITU)) . (TIT)) , and m and orthogonality of Wigner 3-j symbols Q, we get that in the new basis the 
Fourier transform of the multipole matrix is expressed in the following form 

pOO 

U,1' (k) = Ttt {-iy / dRR^jj {kR) g^, {R). (16) 

Jo 

In this way the three-dimensional Fourier transform given by expression ([U]) of a multipole matrix satisfying rotational 
symmetry ([^ is reduced to relevant Hankel transform given by expression (|16|) . 

Calculations performed in this section can be repeated with very minor modification in order to reduce the inverse 
three-dimensional Fourier transform of a multipole matrix to the one-dimensional Hankel transform. We omit the 
derivation giving only the result for the inverse Fourier transform of a multipole matrix satisfying rotational symmetry, 


9i,i' (R) — 


27r2 


dk k'^jj {kR) gf ,, {k). 


0 


(17) 
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V. MULTIPOLE MATRIX IN POSITIONAL SPACE 


To calculate the inverse Fourier transform of 


G(k) 


J Lmd'm' 


with the use of equation (113, gfii {k) is needed. We 


calculate it by means of transformation m and expression (I3|) for Gimi'm' The calculations yields 

gl, (fc) = dTT {2j + 1) [(21 + 1) (21' + ( 0 0 0 ) ' 


(18) 


It is worth noting that the Wigner 3 — j symbol is not vanishing only when the 1,1', j satisfy triangular inequality, 
|1 — I'l < j < |1 + I'j. In calculations of the above formula we used the following property of the spherical harmonics, 
=^n^.o((2J + l)/4^)^/^ 

With the above expression for gj (k), equation (fT7l) yields 


with 




poo 

9iv {R) = / dk jj (kR) ji (ka) jv (ka) 

Jo 


(2j + 1) [(21 + 1) (21' + 1)]'/^ ( J J 


(19) 


( 20 ) 


It demands to calculate integral of three spherical Bessel functions. That has already been considered in the literature 
0- For R > 2a, the integral is given as follows 


dk ji^ (kR) ji (ka) ji> (ka) 


1-3/2 


/ay+^'+i F (i + 1 + 1') 
\r) r(| + l)F(| + l')^ 


with Euler Gamma function F (x). Whereas for R < 2a we have 

pOO 

/ dk jj (kR) ji (ka) ji> (ka) 


t^/^R 


2a 

— 

a ^ 

7r3/2 

/R 

2a 

[a 

7r3/2 

i?2 


4 F 3 


-I-I' 1 + 1-1' 1-1 + 1' 2 + 1 + 1' 1 3-7 4+7 R^ 


I3j,i,i' 4 -F 3 


2 2 2 2 2 2 

j — I — V — 1 j + I — I' j + I' — I 1 + 1' + j + 1 1 + j j 3 


2 

r +1 +1 


2a o2 2 ' 2 

with coefficients and given by 

r(¥) 


1 i-r , I'-i , i + i' + i 3 ^ J- 
,1 + ^,1 + ^,!+ ^ ’ 2 ’^“ 2 ’^' 


R^ 

J + 1 


2 ’ 2’ 2 4a2 


2 ’ 4a2 


/3,-= 2-3/2 


r(i±^)r(i±^)r(2^)’ 

r(i-j)r(i±i^) 

F (1 + y(i + F (1 + r (I + j)' 


ij 


— _9-7/2_ 


r(i-i) 


r(^)r(^)r(2 + 2±i)' 

Symbol 4+3 stands for hypergeometric function 


( 21 ) 


(22) 


( 23 ) 


4F3 (q?!, a2, a3, a^; / 3 i, ( 32 , Ps', x) 


00 


E 


(ai)k («2)fc (a3)fc ia4)k 

(I3i),(P 2),(P3), kl 
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where ia)j. denotes Pochhammer symbol defined by 

(a)Q = 1 for k = 0 

{a)f. = a (a + 1 )... (a + fc — 1) for k = 1,2,... ' 

For nonoverlapping configurations R > 2a, the multipole matrix elements are given by formulae (|19|l . (j20l) . and 
dm). Therefore in this regime gj^, (R) is proportional to 

For overlapping configurations, R < 2a, the multipole matrix elements are given by formulae (IT^ . ((20l) . and (|2^ 
with equations (1231) . For small I and I', analysis of the poles in Euler gamma functions and in Pochhammer symbols 
(which may also be expressed by Euler gamma functions) reveals that for overlapping configurations g^, (R) is a 
polynomial with respect to R. We observe that the degree of the polynomial is ^ + F + 1. For example 


A (/i < 2a) = 

(24) 

, , , (-'ia'^R + R^f 

(25) 


VI. SUMMARY 

In this article we calculate the multipole matrix elements of Green function for Laplace equation. The elements 
are defined by equation ©• The expression for nonoverlapping configurations, i.e. for R > 2a, is known in the 
literature and can be inferred e.g. from the reference Q- The new contribution is the calculation of expression ([T]) for 
overlappiM configurations, i.e. for R < 2a. In this case one can find related considerations for the lowest multipole 
numbers [8|. 

It is worth mentioning that the method which we use to calculate multipole matrix for Laplace equation can be 
generalized for other cases, because the considerations rely on simplification of the Fourier transform. For example, the 
method can be generalized to the multipole matrix elements of Green function for Stokes equations in hydrodynamics 
[IMl. In this case overlapping configurations of multipole Green function for the lowest multipoles have been 
considered recently in the literature 0 - 
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